Elastic scattering losses from colliding BEC's 
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Bragg diffraction divides a Bose-Einstein condensate into two overlapping components, moving 
with respect to each other with high momentum. Elastic collisions between atoms from distinct 
wave packets can significantly deplete the condensate. Recently Ziri et al. (Phys. Rev. Lett. 94, 
200401 (2005)) have introduced a model of two counter-propagating atomic Gaussian wavepackets 
incorporating dynamics of the incoherent scattering processes. Here we study the properties of this 
model in detail, including the nature of the transition from spontaneous to stimulated scattering. 
Within the first order approximation we derive analytical expressions for the density matrix and 
anomalous density which provides excellent insight into correlation properties of scattered atoms. 
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I. INTRODUCTION 

A light-induced potential applied to a Bose-Einstein 
condensate can be used to make high momentum daugh- 
ter BEC wavepackets, which propagate through the par- 
ent condensate 0,0, ~ m ■ Such technique have been used to 
4 , to study the coherence properties 
6 1 , and to study nonlinear four- wave 
mixing (4WM) of coherent matter waves S H Ell- 
in the process of a collision between two condensates, in- 
evitably, some atoms would scatter away, forming often a 
noticeable hallo around the region of collision. The hallo 
at the level above noise indicates a profuse elastic scatter- 
ing losses, and has been observed in various experiments, 
for example [HUEl. 

A scheme of such a process is illustrated in Fig. ^ 
where the collision of two condensates having mean 
wavevector ±Q is shown in the momentum space. The 
condensates are denoted by large dots and marked 
as ip±Q- Atoms from these two counter-propagating 
wavepackets undergo elastic collisions and can be scat- 
tered out into all the modes permitted by the energy and 
momentum conservation (a pair of such states is marked 
with small dots). These modes lie within the three di- 
mensional shell, which in Fig.^is represented by the gray 
ring. 

The elastic scattering loss is not accounted for by the 
Gross-Pitaevskii Equation (GPE). Recently, some modi- 
fications of the GPE were proposed in order to incorpo- 
rate this loss into the mean field dynamics. This has been 
done either within the slowly varying envelope approxi- 
mation by adding an imaginary part to the scattering 
length 13], or by including a stochastic component to 
mimic quantum field fluctuations [hH Il5| . Interesting re- 
sults were also obtained using method based on the field 
theory formulation in the lowest order of the perturba- 
tion approximation 

[lflllZl- So as to test the validity of 
various approximate schemes an exact (nonperturbative) 
model of collisional losses was proposed 18] . It is valid 
both in the regime of spontaneous initiation and Bose en- 
hancement. This model assumes spherical non-spreading 
Gaussians for the colliding wave packets and is capable 
of treating the number of scattered atoms as well as their 



statistical properties through the higher order correlation 
functions. 

Here we study this model in more detail, including the 
nature of the transition from spontaneous to stimulated 
scattering and statistical properties of scattered atoms. 
We also provide a link to previous results obtained within 
the perturbation approximation. 




FIG. 1: Collision of condensates in momentum space. Big 
dots denoted as ip±Q represent colliding condensates. Small 
dots represent scattered atoms. Dashed lines represent ani- 
hilation of atoms, solid lines represent creation of atoms; it 
should help to imagine anihilation of two atoms, one from each 
of the condensates, and creation of two atoms in such states, 
that the momentum and energy of the system is conserved. 



II. THE MODEL 

A system of Bosons interacting via contact potential is 
described by the Hamiltonian 

f fi 2 V 2 
H = - d 3 rW(r,t)— — *(r,t) 
J 2m 

+ | J d 3 r&(r,t)¥(r,t)*(r,t)*(r,t), (1) 
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where ^(r,t) is a field operator satisfying equal time 
bosonic commutation relations, m is the atomic mass 
and g determines the strength of the inter-atomic inter- 
actions. Since the Hamiltonian is of the fourth order 
in the Heisenberg equation governing the evolution of 
the field, 



2yj2 



h'X7 
2m 



•#(r, t) + g¥(r, t)$(r, i)*(r, t), (2) 



is nonlinear and thus, in general, analytically and nu- 
merically untractable. However, for some physical sys- 
tems, a Bogoliubov approximation can be applied lead- 
ing to linear Heisenberg equations. The idea underly- 
ing this approximation states that for some cases the 
field operator might be split into two parts: ip and 5. 
The first contribution describes macroscopically occupied 
modes, where the fluctuations are usually small; hence 
its operator character might be dropped (i/j becomes a c- 
number wave- function satisfying GPE) . The second part 
5, called, "above mean field part", representing fluctua- 
tions, requires full quantum mechanical treatment. Since 
the initial state of the system consists of two counter- 
propagating atomic wave-packets and the "sea" of unoc- 
cupied modes, Bogoliubov approximation can be applied. 
The splitting of the bosonic field is performed in the fol- 
lowing manner: 



tf(r,t) = iP Q (r,t) + i>- Q (r,t) + 5(r,t), 



(3) 



where the subscript ±Q denotes the mean momentum 
of the colliding condensates and iPq(t, t) + ^_Q(r, t) sat- 
isfies the time-dependent GPE. Upon inserting Eq. © 
into the Heisenberg equation one obtains on the right 
hand side several terms, of which, in the spirit of the 
Bogoliubov approximation, we keep only those up to the 
first order in 5, 



iUd t S(r, t) 

h 2 v 2 



(4) 



2m 



2g\^ Q (r,t)\ 2 + 2g\^ Q (r,t)\ 2 



+ 2gi(;* Q (r, t)^_ Q (r, t) + 2gi>*_ Q {r, t)^ Q (r, i)J <5(r, t) 

+g [2^j Q (r, t)iP- Q (r, t) + ^ (r, t)+^_ Q {v, t)} # (r , t) . 

The above equation can be simplified using follow- 
ing arguments. First, we expect that for as long 
as the mean kinetic energy of the scattered atoms 
(h 2 Q 2 /(2m)) is much larger than the interaction en- 
ergy (gn, where n is the mean density of the atoms 
in the condensate), the mean field energy terms 
[2 ff |V Q (r,i)| 2 + 2 5 |^- Q (r,i)| 2 + 2^(r,i)V-Q(r,t) + 
2gtp*_q{Y, i)tpQ(r, t)], can be dropped, leading to 



ihd t S(r, t) = — 



h 2 V- 

2m 



(5) 



To make further simplifications more transparent, we 
visualize effects of various terms in the momen- 
tum space using schematic picture. The first term, 
2g6i(r ) t)ipQ(r,t)ilJ-Q(r,t), is a source term. It corre- 
sponds to annihilation of two particles from counter- 
propagating condensates and creation of two particles in 
the field of the "sea" of scattered atoms. An example 
of such a process is illustrated in Fig. ^ The second 
term, g$(r, t) [^g(r, t) + iP^_q(t, t)] , as shown in Fig. 
describes annihilation of a pair of particles from the same 
condensates (either ipQ or iP~q) and creation of two par- 
ticles in the 8 field. As seen from Fig. [3 this process 
is non-resonant, leading to violation of the conservation 
laws. 




FIG. 2: Visualisations of ^(r.i) [^(r.f) + ip- Q (r, i)]term. 
Figure a) represents the term proportional to ip-Q- Two 
atoms are anihilated from ip-Q and scattered to such states 
that the momentum is conserved. Analogeous plot b) is a 
visualisation of term proportional to ipg. 



The above arguments imply that the Bogoliubov equa- 
tion JSJ can be simplified by dropping the off-resonant 
term and leaving only the source term. The evolution 
equation for the 8 field finally becomes 



ihdtS(r, t) 



h 2 V 

2m 



+g [2^ Q (r, t)i>- Q (r,t)+ ^ (r, t) + (r, t)\ tft ( r , t ). 



S(r, t) + 2 5 V Q (r, i)^-g(r, t). 

(6) 

To simplify the dynamics even further we define and 
compare three characteristic timescales of the problem. 
Let a and N/2 be the width and number of particles 
in each of the colliding condensates. Then, the col- 
lisional time, (the time it takes for each wave-packet 
to pass through its colliding partner) is defined by: 
tc — (ma) /(HQ). Another characteristic timescale of 
the problem is the linear dispersion time (time of the 
spread of the wave-packet due to kinetic energy term), 
tLD = ma 2 /h 9J. In the same manner the nonlinear 
dispersion time, which is the time of ballistic expansion 
in Thomas Fermi approximation |19| . can be defined: 
t^D = v Tr 3 / 2 ma 5 /gN. The dynamics of the system 
depends on the relations between these timescales. It is 
convenient to introduce the following dimensionless pa- 
rameters: t LD /t c — @ and (t LD jt ND ) 2 — a. When the 
number of elastically scattered atoms is small in compar- 



3 



ison with the total number of atoms in both wave-packets 
and both linear and nonlinear dispersion timescales are 
much longer than the collision time {{thD/tc) = /9 > 1 
and (tND/tc) 2 — P 2 / a 2> 1), we can neglect the change 
in population and shape of the macroscopically occupied 
functions iP±q(t, t) during the collision. In our model we 
assume Gaussian shape of the colliding condensates 



exp 



exp 




ihtQ 2 



2m 



(7) 



where r = (x%, X2, £3). Although the assumption that 
the Gaussians have spherical symmetry is not the most 
general, it makes the problem numerically tractable and 
allows to get insight into all relevant physical aspects of 
the system. 

Upon substitution (J7J into JBJ, and then rescaling the 
variables r/er — > r we obtain: 



0d t S(r,t) 



1 



A6(r,t) 



ae 



ipt P(r,t). (8) 



Next we switch to dimensionless field operators substi- 
tuting 5(r, t)a 3 ^ 2 — > 8(r,t). Notice, that assumption of 
spherical symmetry of the colliding wave-packets imposes 
such symmetry on Eq. JHJ. Thus it is convenient to de- 
compose 6 into the basis of spherical harmonics 



S(r,t) = R n,l(r)Ylm( 



ln,l,m(t), (9) 



where a, n j tTn are annihilation operators for a particle in 
the mode described by n,l,m quantum numbers. There 
is still a freedom of choice with regards to the set of 
orthogonal functions R n .i(r). As we shall see below a 
good candidate is a set of eigenfunctions of spherically 
symmetric harmonic oscillator, 



2n\a ~ 3 



c 



(10) 

where L"n 2 (x) is the associated Laguerre polynomial [2(j 
and ao, a harmonic oscillator length, is an auxiliary free 
parameter that can be chosen to minimize the computa- 
tional effort. Notice, that the choice of orthogonal basis 
preserves bosonic commutation relations for annihilation 
and creation operators a n i tm . 



i a n,l- 



Sw 



11 Oifi'l' 



0, 



[ a n,i,m> a n'l' 



0. 



(11) 

(12) 



Substituting Eq. JJJ into Eq. (JSJ) and making use of the 
completeness of the basis functions we obtain the equa- 
tion governing the evolution of operator a ni ; l7Tl , 



id+a 



0>n,l,m + D n ia n — 1. 



- t2 E 



C n , n ',l%> 1 - m , (13) 



where E n l = (2n + I + 3/2)/ag, B n<l = (E n<l - (3 2 )/2(3, 
D n j = y/n(n + l + 1/2)/ (2/3a%), and 



C n , n ',i = / r dr Rn t i(r) exp(-r )Rnn(r) 
Jo 



r(n + i + |)r(n' + z + §) t_ 
^ r(z + |) 2 r(n + i)r(n' + i) 1 +a ° } 



1 



n+n 



F{-n,-n',l+^,l/at\. (14) 



Here F(a,b,c,x) is a hypergeometric function |2^. No- 
tice that all coupling coefficients are calculated analyti- 
cally and the a n ,i,m operators for different I and m are 
decoupled. Moreover, equations l|13|) do not depend on 
quantum number m. With all these simplifications the 
linear system of equations (jl3|l can be solved numerically. 
The general three-dimensional problem, due to spherical 
symmetry simplifies to the set of one-dimensional ordi- 
nary differential equations. In numerical applications, for 
every / a basis of approximately 64 modes associated with 
quantum number n is sufficient. Moreover, simulations 
show that 50 modes associated with quantum number I 
is usually enough. Thus the whole quantum model can 
be solved numerically on an ordinary PC within one hour 
of calculations! 

We define a vector of operators: 



&i,l,r 



Vj, m (*) = 



\ 



"0,1, 



"1,1, 



and rewrite the equation l|13|) in the compact form: 



i-n v l,m{t) = Aj(i)vj, m (f). 
at 



(15) 



The matrix A;(i) is plotted schematically in Fig. EI 

The solution vj jm (i) can be expressed in terms of the 
time evolution operator: 

(f)=U,(t)vi, m (0). (16) 

The evolution operator satisfies equation 

dU,(i) 



df 



-iA,(f)U,(t) 



(17) 



with the initial condition U;(0) = 1. In our calculations 
we evaluate A;(t) and find evolution matrix Uj(t) by 



4 



\ "'" D 

n+1 1 

B , \ 

n,l 




^ t! c„.,„ 


-- e t2 r 


-D 

n+1,1 

\ ~B , \ 

n,l 



FIG. 3: Visualisation of A;(t) matrix from Eg. 1151 . It is 
tri-diagonal in coupling a ni ; im 's with a n / ; , m 's and a\ ll _ m 's 

with aL , _ m 's. The coupling between a n ,i,m's and a , , _ m 's 
— t 2 

is determined by — §e C n , n ',i- In general, all s are 

non-zero, which is represented by grey area. 



solving Eq. I|17f) . The matrix U;(t) uniquely determines 
full operator dynamics of <5(r, t) and hence is sufficient to 
obtain any expectation value of the variables referring to 
the system. 

Equation JSJl can be treated alternatively in the Fourier 
domain. Upon defining the Fourier transform in the form 



6(r,t) 



d 3 k exp (i(3kr - i(3k 2 t/2) <5(k, t), 



we can rewrite Eq. JHJ as 



id t S(k,t) 



8tt 3 / 2 



e / d k exp 



x exp ( —i(3t 1 - 



k 2 + k' 



(18) 



We took this particular form of the Fourier transform 
(with factor f3) so that the wavevector of the moving 
condensate is equal to unity. 



III. FIRST ORDER PERTURBATION THEORY 

The advantage of using Fourier transform lies in the 
fact, that the equation (|18l) can be used in a natural 
way as a basis of the perturbative expansion. Assum- 
ing a general solution of S(k, t) in a perturbative form, 
<5(k, t) — (k, t) + 5^> (k, t) + . . ., we derive a recurrent 
relation, which in the lowest order gives 



/? 2 (k + k' 



i a t ^)(k,<) = ^e- t2 |d 3 fc'exp 
x exp f-if3t ( 1 - m(k',t). (19) 



The first order equation (|19fl can be integrated and con- 
sequently we can evaluate observable quantities like num- 
ber of scattered atoms or the correlation functions ana- 
lytically. 



A. Density matrix 

In this subsection we find the approximate expression 
for the density matrix of scattered atoms, 



p(k 1 ,k a ,t) = <5t(k 1 ,t)*(k 2 ,t)), 



(20) 



within the first order perturbation theory. We perform 
the calculations in the momentum space. 

Our system is initially in the vacuum state |f2), that is 
<j(k,0)|fi) = 0. Hence, due to the commutation relation, 
[<5(k, t), <5t(k', t)] — S 3 (k — k'), at initial time we have 



(<5(k,0)^(k',0)) = <5 3 (k-k'). 



(21) 



Then, upon substituting Eq. (|T$)l into (|2U|l . and using 
Eq. we obtain 

P(^M = ^ fdr Tdr'exp^-r' 2 ) 
J o Jo 



x exp 



64tt 3 
P 2 (k\ + kl) 



if3(Tkl-T / k 2 )/2 



d k exp I — 



/3 2 (fc 2 +k(k!+k 2 )) 



xexp(i/3(r-r')(l-/j 2 / 2 )) 



(22) 



We show how to handle the above integral in Appendix 
A. For (3 1 and (3\k-i +k 2 | ~> 1 the density matrix 
reduces to 



a 2 /3V2n 
32tt 2 



p(ki,k 2 ,i) 
x J dr J dr' exp ^— 



exp 



/5 2 |ki 



(r-r') 2 |k 1+ k 2 | 



x exp ( — t 2 - r' 2 + i/3r ( 1 - 
x exp I —ifir' I 1 — 



|k!+k 2 | 2 



kj 
2 



ki+k 2 | 2 



A" 2 



(23) 



This expression, although it still contains twofold inte- 
gral over time, is simple enough to give the distribution 
of scattered atoms, and provide a further insight into co- 
herence properties and dynamics of scattering process. 



B. Final momentum distribution and coherence 
properties of scattered atoms 

From the last expression in the previous section 
(Eq. H23|) ~) the density matrix becomes negligible if 



FIG. 4: Configuration of vectors ki and k2 as in case a). The 
vectors are of the same length, with angle <f> between them. 
The dashed circle corresponds to k = 1 and the grey area is 
the scattering shell. The length of each of the vectors differs 
from unity by Afc. 



FIG. 5: Configuration of vectors ki and k2 as in case b). The 
vectors are parallel. The dashed circle corresponds to k = 1 
and the grey area is the scattering shell. The lengths of ki 
and k2 differ from unity by ±Afc. 



any of the arguments /3|ki — k^l, (3 yt 
or (3 ( 1 — l kl + k2 l ^ ) , appearing in the exponential 



|k!+k 2 | 2 _ ki 

8 2 



8 2 

functions, are large (remember that we work under con- 
stant assumption (3 3> 1). Hence the only regions where 
we expect nonvanishing value of p(ki,k2,i) are those, 
for which wavevectors ki and fc 2 are close to each other 
(|ki — ki2 1 <C 1), and they both have the length close to 
unity (notice that in our units this coincides with the 
length of the wavevector of the colliding condensates) . It 
is interesting to examine two cases: 

a) Equal length case; fci = fc 2 = 1 + Afc. In this case 
(see Fig. the vectors ki and k 2 are of equal length 
and form an angle <f>. The properties of the density ma- 
trix strongly depend on the distance from the resonance 
surface Afc. For small values of Afc and small angles <f> 
the density matrix can be approximated by the Gaussian 
function 



p(Afc, (b) = .A • exp - 



(3 2 <p 2 (3 2 Ak 2 



The constant A = 



a 2 p 
128^ 



1 + ^ArcTan 



( i 

\2\/2 



(24) 



was 



found evaluating (|23[) analytically on resonance (ki = k 2 
and |ki| = 1). Notice that Eq. for <f> = repre- 
sents the density of the cloud of scattered atoms. The 
cloud has spherical symmetry, following the symmetry of 
Eq. (JSJ). It takes a form of a shell localized around k = 1 
(which corresponds to wave- vector Q in our units). Its 
width is of order to 1/(3. 

b) Parallel vectors case; ki || k2 (see Fig|SJ). Wavevec- 
tors are shifted away from the resonance in the opposite 
direction (ki = e(l — Afc) and k2 = e(l + Afc)). In this 
case we also find a Gaussian fit 



p(Afc) = A ■ exp 



--(3 2 Ak 2 +i(3Ak 



This concludes our analysis of the density matrix. Now 
we move on to the first order correlation function 



3i(ki,k 2 ,i) = 



(Jt(k 1 ,t)<S(ki,t))((Jt(k 2 ,t)J(k 2) *)) 



We find that in the case a) 



5i(Afc» = exp 



while in the case b) 



(3 2 (b 2 



(26) 



0i (Afc) = exp 



--(3 2 Ak 2 + i(3Ak 
6 



(27) 



(25) 



By comparing Eq. 127|) describing the coherence in the 
radial direction - <?i(Afc), and Eq. (|24|l describing the 
density p(Ak, <f> — 0) it seems reasonable to assume that 
the radial coherence length of the shell of scattered atoms 
is of the same order as its width. Following this obser- 
vation we present a very rough estimate of the number 
of atoms necessary for the bosonic stimulation to occur. 
First we introduce the concept of the coherence volume 
associated with each scattered atom. It is determined by 
the first order correlation function. If we fix ki and vary 
k 2 we get a well defined peaked structure around k 2 = ki . 
The size of this structure gives the coherence volume as- 
sociated with the single atom (see Fig. As we have 
already shown, the radial extend of this structure is of 
order of 1/(3. From Eq. (|26|l we deduce that the angular 
coherence length is of the same order 1/(3. Hence we esti- 
mate the coherence volume to be equal to 1//3 3 . Next we 
remind ourselves that the volume V in the momentum 
space accessible by the scattered atoms is spherical shell 
determined by the diagonal part of the density matrix. 
In our case V — An/ ' (3. When the product of the number 
of scattered atoms times their coherence volume is com- 
parable with V, the coherence volumes associated with 
distinct atoms start to overlap (see Fig. EJ. Then, the 
wave nature of matter comes into play, and interference 
effects would enhance the scattering - bosonic stimula- 
tion |2l| . Using the above estimates, the critical number 
of atoms is equal to 4tt(3 2 . In the next section, we derive 
more practical criterion for bosonic stimulation to show 
up in the collision of the condensates. 



6 




FIG. 6: The artistic representation of the atomic coherence. 
Final states of scattered atoms lie on the spherical shell and 
the coherence volume of each atom is marked as a dark region 
in this shell. Case a) corresponds to weak scattering and in 
b) individual coherence volumes start to overlap. 

C. Number of scattered atoms for t^$> 1//3 

The most straightforward observable quantity, the 
number of scattered atoms as a function of time can be 
expressed in terms of the trace of the density matrix 

S(t) = J d 3 k p(k,k,t). (28) 

Using Eq. we get: 

5 (i)( t ) = _^£ = r k 2 dk f At f At' exp (-t 2 - t' 2 ) 
4V2^ Jo Jo Jo V 1 

x exp (i/3(l - k 2 )(r - r') - y (r - r') 2 ) . (29) 

Unless k is close to unity the integral over t and r' will 
vanish due to the rapidly oscillating phase factor 13(1 — 
k 2 )(r — t'). Thus it is convenient to change the variables 
k = 1 + § and expand the integrand in the lowest order 
of x obtaining 

s(1) = aV2^ f°° dx f* dr f*^ / 2_ T /2N 
8tt J_oo Jo Jo 

x exp^i2a;(T-T')-i(T-r') 2 j • (30) 

Next, evaluating integral over x we get 

S (D {t) = ^1 f dT [\t'8(t-t') 
8 Jo Jo 

x exp(^-r 2 -r' 2 -i(r-r') 2 ), (31) 
and finally the number of scattered atoms is equal to 

S (1) (i) = ^rerfV2). (32) 



In Appendix B the same result is obtained within a clas- 
sical model of colliding hard spheres with cross-section 
equal to 87ra 2 . 

An important remark concerns the condition for 
bosonic stimulation. As we mentioned in the previous 
section it occurs when S > Air/3 2 . By comparing this 
result with Eq. (with t — > oo), we get a//3 > 2. 

This is a simple condition, which translated into experi- 
mental parameters provides the criterion for bosonic en- 
hancement in the collision of two condensates. To verify 
this condition we evaluate the total number of scattered 
atoms up to the third order 

5 <.» +5 <*, +5( 3, = 4( 1 + c ,5! + ^), (33) 

where ci and C2 are numerical coefficients. This suggests 
the functional dependence of the form S — ^j^f ffk)) 
where /(0) = 1. Our previous numerical results 
|18| , obtained both in perturbative and non-perturbative 
regimes confirm this suggestion. In Fig0we plot the log- 
arithm of / versus (a//3) 2 and obtain a straight line with 
high accuracy. Hence we see that indeed the bosonic 
stimulation (manifested by the value of f(a 2 /f3 2 ) be- 
ing larger than unity) occurs for a/j3 ~ 2. More- 
over, function / has an exponential form f[(ce/f3) 2 ] — 
exp[0.085 (a/P) 2 ]. Using above formulas one can both 
easily estimate number of scattered atoms and predict 
whether bosonic enhancement would be present in the 
particular physical realization. 




a 2 /p 2 

FIG. 7: Elastic collision loss as a function of a 2 /f3 2 . The 

total number of scattered atoms is given by 5 = Trfg / ^fh^ • 

The figure shows that, ln(/) is a linear function of a 2 / f3 2 with 
slope equal to 0.085. 
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D. Number of scattered atoms for i< 1, 
Characteristic timescales of the signal build up. 

In order to analyze the dynamics of (S'S) on a very 
short timescale, we consider the diagonal part of the den- 
sity matrix 



p(k, k, t) = jf * dr jf* dr' exp (-r 2 - r' 2 ) 



x exp 



^_(W) 2 fc 2 +i/3(r __ r , )(i _ fc2) ^ 



(34) 



Notice that for very short times (t <C 1), all the terms 
proportional to squares of r and r' in the exponents ap- 
pearing in Eq. can be dropped. With this simpli- 
fication the time integrals can be evaluated analytically, 
giving 



p(k,M) 



2 /V2^sin 2 (/3t(l - k 2 )/2) 



(35) 



32tt 2 - fc2)/2) 2 

The time evolution of the density of scattered atoms in 

a: -2-^ _ hi -2 h 

□ □ 

-2-10 1 2 ,-2 -10 12 



c)- 



-1 



o 



FIG. 8: The density of atoms in momentum space as given 
by Eq. for a = 10 and j3 = 20. The halo of scattered 
atoms shrinks, forming a shell around k = 1. The four frames 
correspond to increasing times a) t = 0.07t c , b) t = 0.1i c , c) 
t = 0.15i c and b) t = 0.4t c , where i c is the collisional time. 

k x -k y plane is pictured in FigEl It shows that atoms, 
which initially scatter isotropically in the momentum 
space, eventually form a well pronounced halo around 
the k = l. 

To relate Eq. ffify to the number of scattered atoms, 
we need to perform the integration over k. Evaluating 
angular integral explicitly and changing the variables x = 
(3t(k 2 - l)/2, we get 



S{t) 



a 2 V^ f 1 



8- 



sin 2 x Fix ~ 

pt/2 x 2 y pt 



a) /3t <C 1; we set the lower limit in the integral equal 
to zero and use the approximation y|jf 
obtain 



1 



'22 to 
/Si lo 



S(t) 



/2tt 



-U 



dx 



sin 2 x 



a 



8tt " V fit Jo x 2 v 4v^ 

b) /3t 3> 1; we set the lower limit in the integral equal to 



-co and use the approximation \/M 



1 ~ 1 to obtain 



S(t) 



2tt 



8tt 



da; 



• 2 2 

sin x cr 



2tt 



Let's discuss the validity of the above result. Notice that 
we introduced approximations describing two body inter- 
actions by a contact potentials (renormalized delta func- 
tion potential) . On the other hand any realistic potential 
should have a natural ultraviolet cut-off k c . If we extend 
the upper limit only up to k c (instead of extending it to 
oo) we find in the case a) discussed above a narrow win- 
dow of very short times when the number of scattered 
atoms grows quadratically in time. On later times dy- 
namics will not be affected by the cut-off. Consequently, 
the realistic estimate would predict that the number of 
scattered atoms should be quadratic function of time for 
initial extremely short period, next there is a window of 
time when it behaves like a square root, and finally it 
turns into linear, semiclassical, regime. 



E. Higher order correlation functions 

A lot of interesting information about the quantum 
systems might be obtained from its correlation functions. 
In this section we present and discuss some of their prop- 
erties. We focus on equal-time correlation functions in 
momentum space. This choice is justified by the fact 
that in most experiments performed with Bosc Einstein 
condensates, the common imaging technique is time-of- 
flight absorption, which directly registers the momentum 
distribution of the atomic cloud. 

The normalized n-th order density correlation function 
(n > 2) is defined as: 



5n(kl 



> k n , t) 



(tf(k 1 ,t)...tf(k n ,t)8(k n ,t)...6(k 1 ,t)) 

(SHk 1 ,t)S('k 1 ,t))...(s^k n ,t)S(k n ,t)) ' 



(36) 



We make use of fact that the evolution of 5 is linear, 
which means that 5(k, t) is a linear combination of oper- 
ators <5(k',0) and ^(k^O); hence 

5{k,t)= J dk! (u {k,k! ,t)5{k! ,Q) + V{k 1 k! ,t)P {k! ,Q) 



This formula can be further simplified in the following For such a system Wick's theorem applies. Since the 



two regimes: 



averages are calculated in the vacuum state |f2), where 
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S(k, 0)|fi) = for all k, all correlation functions decom- 
pose into combination of products of a density matrix 
p(kj,kj,t) and anomalous density m(ki, k m , t), defined 

as 



takes the form 



m(k ll k ai t) = ^(ki,t)^(k a ,*)). 



(37) 



This leads to the growing complexity of the correlation 
functions, since in the n-th order, the number of terms 
contributing to the g n is of the order of (2n — 1)!!. For 
example, the second order correlation function 



g 2 (k ly k 2 ,t) 



(fr(k u t)8Hk 2 ,t)6(k 2 ,t)8(k u t)} 



g 2 (ki,k 2 ,t) = 1 



|p(k 1 ,k 2 ,t)| 2 + |m(k 1 ,k 2 ,t)| ; 
/ o(ki,ki,i)p(k 2 ,k 2 ,i) 



(38) 



In the first order, we find an analytical expression for the 
anomalous density 



m(ki,k 2 ,i) 



a/3 2 



■ exp 



/3 2 (k! + k 2 ) 2 



.a/3 2 

-i exp 

16tt v 



/3 2 (k 1+ k 2 ) 2 /3 2 



dr exp ( -t 2 - i0r ( 1 - ^±M 



1 



k\ 



(39) 



erf ( t + -0 ( 1 



The above function reaches its maximum when ki = — k 2 
and both vectors are of unit length (fci = k 2 — 1) due to 
the energy and momentum conservation. 







i i 
opposite momenta 






correlation^^^^^^ 




atom 






/bunching 








i 


Ill 



30 60 90 120 150 

<l>[deg] 



FIG. 9: The second order correlation function for ki = fe = 1 
as a function </> - relative angle between the vectors and a = 
20, (3 = 30. The strong maximum around <f> = corresponds 
to atom bunching. The maximum around <f> = 180° is an effect 
of correlations of pairs of atoms with opposite momenta. 

Combining Eqs (|38|l and 139|) with the expressions for 
the density matrix, discussed before, we can find ana- 
lytic expression for the second order correlation function 
in the first order perturbation theory. In Fig|5]we show 



the second order correlation function for k\ = k 2 = 1 as 
a function of relative angle </> between the vectors. The 
maximum around (f> — corresponds to the detection of 
two scattered particles close to each other. It is greater 
than one due to the bosonic statistics of scattered atoms. 
This effect is called atom bunching. Since the atoms scat- 
ter in pairs of opposite momenta, the g 2 correlation func- 
tion reveals strong correlation around tfi — 180°. It is ex- 
clusively an anomalous density part of the second order 
correlation function that contributes to this maximum. 



IV. CONCLUSIONS 

We studied the nature of the transition from spon- 
taneous to stimulated scattering regime in the collision 
of two Bose Einstein condensates. Within the first or- 
der approximation we derived analytical expressions for 
the density matrix and anomalous density. This enabled 
us to derive the condition for bosonic stimulation effect. 
The characteristic timescales in the dynamics of scattered 
atoms were identified. Finally, we showed that scattered 
atoms feature both bunching as well as the opposite mo- 
menta correlations. 
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VI. APPENDIX A / /3 2 fc 2 , , N A fc 2 

x exp I — + »/3(t - r ) I 1 - — 

In this appendix we simplify the expression for density 
matrix (|22fl upon evaluating the integral over d 3 k in the 
limit /3 > 1. The angular part of this integral can be 
performed analytically 

/ B 2 (k 2 + k(ki + k 2 ))\ Next ' using the idcntit y 2 ( fc i + k l) = l k i+ k 2 

d 3 k exp ( I |ki — k 2 | 2 , we get 



a 2 f3 2 



p(ki,k 2 ,i) 



1 



16tt 2 |ki + k 2 



exp 



/3 2 |kx 



J kdkexp(-i^(T-T')k' 



exp 



dr / dr' exp (— r 2 
- exp - 



^'2 



/? 2 



i/3(Tfc 2 -T'fc 2 )/2 + i/3(r-r')) 



(40) 



The most significant contribution to the above in- 
tegral comes from the region, where the condition 

/3|ki+k 2 | ^> 1 is satisfied. In this region, we can 

2^ 



neglect exp 



[ki+kaj 
2 



in comparison with 



J 



a new variable fc 



[ki+kaj 
2 



Additionally we introduce 
^ and get 



f 00 / 3 \ ( B 2 ( 

J kdk exp f — ^(t — r')k 2 j exp I — — |^fc 



|ki+k 2 | 



(41) 



1 



dx (/31k! + k 2 |/2 + x) exp --(r - r') /?|k x + k 2 | 2 + 4x|k x + k 2 | + — exp (-x 2 /2) 

-/3|k!+k 2 |/2 V O V P // 



4^ 



The condition /3|ki+k 2 | 3> 1, together with the 
presence of the exponential factor exp (— x 2 /2), allow 
to extend the lower limit in the integral over x in 
Eq. I|41l) to — 00, and neglect x in comparison with 
/3|ki +k 2 |/2. Due to the presence of the exponential 



factor exp (— r 2 — r' 2 ) in the equation l|40|) we can make 
another approximation, namely neglect the phase factor 
*(t — r')#g. Having all above approximation applied we 
can evaluate simple Gaussian integral to get finally 



p(ki,k 2 ,i) 



a 2 BV2^ 
32tt 2 



■ exp 



/3 2 |ki-k 2 | 2 

8 



dr / dr' exp 
i] Jo 

kx+k 2 | 2 k 2 



^2 
-T — T 



(T-r0 2 |k 1+ k 2 | 

8 



(42) 



VII. APPENDIX B 



model is meant as a classical counterpart of the collision 



In this appendix we calculate the collisional loss for two 
counter-propagating clouds of classical particles. This 
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described in the main body of the paper. To match the 
conditions used in the paper we assume a Gaussian den- 
sity profile for each of the clouds and we restrict ourself 
to the dilute gas limit (naoa < 1), where oo is a cross- 
section for a collision of two particles and a is the radius 
of the cloud. In this limit we neglect secondary collisions. 
We also ignore the depletion of the clouds; this approxi- 
mation is valid as long as the fraction of scattered atoms 
is small. When the collision is described in the reference 
frame associated with one of the clouds, the density of 
this cloud is equal to 



Here, the prefactor "2" on the right-hand side accounts 
for the fact that in every collision two particles are scat- 
tered. For two identical bosons, the cross-section <tq at 
low-energy limit is equal to 87ra 2 , where a is the s-wave 
scattering length. Using the identity 



exp 



erf 



x + 2vt 



dx = \pHa ■ erf 



y/2vt 



m(r) 



N 



2tt 3 / 2 ct 3 



exp 



while the other cloud propagates along the x\ axis with 
velocity 2v, 



n 2 (r, t) 



N 



2 7 T 3 / 2 f 7 3 



exp 



(xi + 2vt) 2 + xl+x\ 



Here both clouds contain N/2 particles. 

The number of particles scattered up to the time t 
equals 



S{t) = 2 x 



d 3 r m(r) / dt' (n 2 (r,t') • 2va ) 
Jo 



and rescaling time (t/tc — > t), where the collisional time 
tc = cr/v, we obtain: 



S(t) 



N 2 a 



2„2 



rf(i\/2) 
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a ir 

IT 



tf(tV2). 



This result reproduces number of scattered atoms ob- 
tained in the quantum model within the first order per- 
turbation theory. 



[1] Yu. B. Ovchinnikov et al, Phys. Rev. Lett. 83, 284 
(1999). 

[2] M. Kozuma et al, Phys. Rev. Lett. 82, 871 (1999). 

[3] J. Stenger et ai, Phys. Rev. Lett. 82, 4569 (1999). 

[4] E. W. Hagley et al, Science 283, 1706 (1999). 

[5] E.W. Hagley et al, Phys. Rev. Lett. 83, 3112 (1999). 

[6] M. Trippenbach et ai, J. Phys. B 33, 47-54 (2000). 

[7] G. Lenz, P. Meystre, and E. W. Wright, Phys. Rev. 
Lett. 71, 3271 (1993), E. Goldstein, K. Plattner, and 
P. Meystre, Quantum Semiclassic. Opt. 7, 743 (1995), 
E. Goldstein, K. Plattner, and P. Meystre, J. Res. Natl. 
Inst. Stand. Technol. 101, 583 (1996), E. Goldstein and 
P. Meystre, Phys. Rev. A 59, 1509 (1999). 

[8] M. Trippenbach, Y. B. Band, and P. S. Julienne, Optics 
Express 3, 530 (1998), 

[9] M. Trippenbach, Y. B. Band,l and P. S. Julienne, Phys. 
Rev. A 62, 023608 (2000). 
[10] L. Deng et al, Nature 398, 218-220 (1999). 
[11] J. M. Vogels, K. Xu, and W. Ketterle, Phys. Rev. Lett. 

89, 020401 (2002) 
[12] N. Katz, R. Ozeri, E. Rowen, E. Gershnabel, and N. 
Davidson, Phys. Rev. A 70, 033615 (2004); N. Katz, E. 



Rowen, R. Ozeri, and N. Davidson, cond-mat/0505762 
[13] Y. B. Band et ai, Phys. Rev. Lett. 84, 5462 (2000) 
[14] J. Chwedehczuk, M. Trippenbach and K. Rzazewski, J. 

Phys. B, L391, 37 (2004) 
[15] A.A. Norrie, R.J. Ballagh and C.W. Gardiner, Phys. Rev. 

Lett. 94, 040401 (2005), see also A. Sinatra, Y. Castin 

and C. Lobo, Journal of Modern Optics 47, 2629 (2000) 
[16] R. Bach, M. Trippenbach and K. Rzazewski, Phys. Rev. 

A 65 063605 (2002) 
[17] V. A. Yurovsky, Phys. Rev. A 65 033605 (2002) 
[18] P. Ziri et ai, Phys. Rev. Lett. 94, 200401 (2005) 
[19] Y. Castin and R. Dum, Phys. Rev. Lett. 77, 5315-5319 

(1996); Yu. Kagan, E. L. Surkov, and G. V. Shlyapnikov, 

Phys. Rev. A 55, R18-R21 (1997) 
[20] M. Abramovich, I. A. Stegun, "Handbook of Mathemati- 
cal Functions With Formulas, Graphs and Mathematical 

Tables", Dover Publications 1974. 
[21] R. P. Feynman, "Feynman Lectures On Physics", Volume 

3, Addison Wesley Longman (1970). 
[22] L. Mandel and E. Wolf, "Optical Coherence and Quantum 

Optics", Cambridge (1995). 



